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Abstract 

A numerical implementation of the transition state theory (TST) is 
presented which can be used to calculate the attempt frequency /o of arbi- 
trary shaped magnetic nanostructures. The micromagnetic equations are 
discretized using the finite element method. The climbing image nudged 
elastic band method is used to calculate the saddle point configuration, 
which is required for the calculation of fo- Excellent agreement of the 
implemented numerical model and analytical solutions is obtained for sin- 
gle domain particles. The developed method is applied to compare fo for 
single phase and graded media grains of advanced recording media, fo 
is predicted to be comparable if the maximum anisotropy is the same in 
these two media types. 
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1 Introduction 



A detailed knowledge of the thermal stability is of utmost importance for mag- 
netic nanostructures for various applications ranging from hard disc media, 
MRAM devices to permanent magnets. A well established tool which origi- 
nated in chemical rate theory for determining the thermal stability relies on the 
transition state theory (TST). In the TST the thermal stability of a minimum 
energy state Mi is determined by the application of the Arrhenius-Neel law, 
T = ^6*°^^ J where the energy barrier AE separates the minimum energy state 
Ml and the saddle point configuration 5*1. The attempt frequency /q limits the 
probability for reversal. 

The thermal stability, which requires the calculation of the energy barrier 
and the attempt frequency was first calculated for single domain particles in 
the intermediate to high damping limit (IHD) for an external field applied ex- 
actly parallel to the easy axis [2 . For systems where the symmetry is broken 
either by an oblique external field or by an additional anisotropy various works 
are published giving analytical formulas for the attempt frequency for different 
damping limits. For the IHD limit formulas for the thermal stability can be 
found in Ref |19l [T21 121 S|. For the limit of very low damping a formula for 
the thermal stability was derived by Smith et al [T^. In order to extend the 
calculation of the damping limit to all values of the damping constant Coffey et 
al. [5]and Dejardin et al. |6]have shown that the Mel'nikov-Meshkov formalism 
can be applied to magnetic systems. Analytical formulas for the intermediate 
damping regime are given by Garanin et al |S] . 

All the previous mentioned work were restricted to reversal via homogeneous 
rotation. Inhomogenous states were treated by Braun who calculated the ther- 
mally activated reversal in elongated particles [T]. Besides analytical formulas 
for the attempt frequency he found that the thermal stability depends on the 
domain wall energy in these elongated particles. 

The approach of Braun was extended to soft /hard bilayers, where the saddle 
point configuration can be described by a domain wall at the soft /hard interface 
by Loxley et al [TS] . 

The paper is structured as follows. In section II details of the implementa- 
tion of the transition state theory in the finite element package FEMME |22|is 
given. The implementation method is compared to analytical solutions in sec- 
tion III. An improved method for the calculation of fo is presented in section 
IV. Summary and discussion are given in section V. 

2 Numerical Calculation of the Attempt Frequency 

The numerical implementation is based on Kramers transition state theory |13| . 
which Langer transposed and expanded to multidimensional systems in 1967 
|14| . According to Langers approach the attempt frequency determines the prob- 
ability current for the configuration distribution on the energy surface around 
the saddle point. 
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The following steps are required in order to calculate the attempt frequency 
by the means of the transition state theory: 

1. The configuration of the system at the minimum and at the saddle point 
given as the values of all degrees of freedom. 

2. The value and curvature of the Gibbs free energy surface at the minimum 
and at the saddle point with respect to appropriate canonical variables. 

3. The systems' equation of motion in canonical variables along the energy 
surface to describe the dynamics around the saddle point. 

Knowledge of the above properties allows one to calculate the attempt frequency. 
This can be written as: 



/o = (1) 

2.1 Dynamical factor A+ from the Landau-Lifshitz-Gilbert 
equation 

The dynamical factor A-|_ is obtained by solving the noiseless linearized equation 
of motion. In micromagnetics the equation of motion can be written in the form 
of the Landau-Lifshitz-Gilbert equation as, 



-TT^/x h:,, - (/x h:,,) (2) 

where the magnetic polarization J and the effective H^f / field are continuous 
functions of space. In the following, we discretize the continuous equation of 
motion using the finite element method using the box scheme |22| I10|. We can 
write for the effective field on the node point i, 

with Vi being the corresponding volume of spin i and E being the Total Gibbs en- 
ergy. The structure of the Landau-Lifshitz-Gilbert equation leads to a constant 



length of the magnetic polarization 



Js ■ As a consequence the number 



of independent equation is 2N, where N is the number of spins. In order to 
apply the theory developed by Langer the system has to be expressed with its 
canonical variables [13]. For a system with A^-spins the canonical variables on 
node point i are given by |H1 [T7] , 

Pi = ViJsCOs{0i) 

qi = (pi 
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where the polar angle angle and azimuthal angle 9i are defined as, 



Ji ~ Js (sin Oi cos (pi , sin 9i sin (pi , cos 6i) (5) 
The Landau-Lifshitz-Gilbert equation for spin i in canonical variables reads, 
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JsVi{l + a2) smflj 



1 (-JsV,asin^e,^^-JsV,sine,§f\ 



dE 



dOi ~^ sine 



^ s *' i 
a dE 



Let US denote 7721-1 — Pi and 7721 = Qi- We can write the equation of motion 
in the form, 



'2i-l 



dt 



dt 



/2i-l(77l, ...,?72Af) 



(7) 



dt 



dt 



= = ■■ /2»(?7l, ■■■,V2n) 



(8) 



Transition state theory considers the linear dynamics around the saddle 
point. Linearizing the right hand side around the saddle point, we obtain: 





+ 


/ 9/1 


dfi \ 

dri2N 


sp 


diiN 
\ 9vi 


dfiN 
dri2N 1 



sp 



\1l2N — Vsp,2N / 



(9) 



The first derivatives of the energy with respect to the spin components at 
the saddle point is zero, thus the first term on the right hand side is zero as 
well. The second term we define as H^y^'. 



H, 



dyn 



( 9h 



. df2N 

\ dm 



I 



df2N 

dr)2N 



(10) 



The last term we define as vector — Vk — flsp.k- The time derivative of 
equals the left hand side of Equation |9] So Eq. [9] can be written as: 



— W Hdynl^ (11) 

This is a linear system of differential equations which has 2N eigenvalues 
with corresponding eigenvectors i^'^. The solutions are of the form: 

iTk = iTk^e^^' (12) 
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With the exception of one eigenvalue, all other eigenvalues of the matrix 
Hdyn at the saddle point are negative. This single positive eigenvalue is the A+ 
that we need for the calculation of the attempt frequency. 

For the numerical calculation of H^yn, the function f = fi has to be derived 
with respect to the coordinate x = rj^. This was calculated numerically using 
finite differences and a seven point stencil method. 



-f{x - 3h) + 9f{x - 2h) - 4:5f{x - h) + 45f{x + h) - 9f{x + 2h) + f{x + 3h) 



The step size used in the numerical results was h = 0.005 



(13) 



2.2 Statistical factor Qq from the Hessian matrix 

The statistical factor Hq relates the curvature of the total Gibbs energy with 
respect to the canonical variables at the saddle point and at the minimum. The 
curvature is obtained by calculating the second derivative of the energy, 



H 



V 



dhi 

din 



— \ 

N-2 



dh2 



drii 



dri-2 



dh2N-2 

dri2N-2 / 



(14) 



a-'E 

dri2N-2dri 



, a^E 

\ drildV2N- 



dri2N -2dri2N - 

rio is the square root of the ratios of the determinants of the Hessian matrices 
at the minimum and the saddle point: 



/ det Hr, 



\detH, 



(15) 



sp I 



It is important to note, that the theory of transitions as developed by Langer 
is only valid for canonical variables. Suppose that instead of the canonical 
variables (pi, Qi), the polar angles (0^,0^) are used to describe the system (see 
e.g. Ref |TS]) and its derivatives gf-J|-. Then the obtained /o is correct only 
if the magnetic state at the minimum and the saddle point fulfills Oi — 11/2 for 
all spins. This prerequisite is, for example, fulfilled for a single domain particle 
having the easy axis pointing along the x-axis and an external field applied along 
the y-axis as in Ref [18]. The same argument applies for the calculation of A+. 

For the numerical calculation of the derivative the same finite difference 
scheme was used as described in the previous section. 
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2.3 Climbing Image Nudged Elastic Band Method 



The calculation of the statistical prefactor flo and the calculation of A_|- require 
knowledge of the saddle point configuration. A very accurate calculation of 
the saddle point configuration is essential for the calculation of the attempt 
frequency. The nudged elastic band method was used in order to calculate 
the saddle point configuration [201 [7] . It was found that the accuracy of the 
nudged elastic band method is not sufficiently high in order to calculate reliably 
the attempt frequency. In order to improve the accuracy of the saddle point 
determination, a climbing image elastic band method was implemented [llj. 

3 Validation of the method 

3.1 Small magnetic cube 

In order to test the numerical calculation of the attempt frequency, we com- 
pared the analytically obtained attempt frequency for a single domain particle 
with uniaxial anisotropy to the numerically obtained one. For the numerical 
model the following parameters were used: damping constant a = 1 , anisotropy 
constant Ki — l.OM J/m^ , magnetic saturation polarization Js = O.ST and ex- 
change constant A — lOpJ/m. The model size is 0.6 x 0.6 x 0.6 nm''. 12 finite 
elements were used to discretize the cube. The easy axis of the cube was as- 
sumed to be parallel to the x-axis. The external field was applied parallel to the 
y-axis. Results were compared to those obtained using the analytically derived 
expression for the attempt frequency given in Ref. [TF. 

In Fig. [l] (B) two magnetization states of the homogeneous cube are shown: 
on the left the initial state is depicted, with magnetization pointing close to 
the the easy axis. On the right the magnetization is shown at the saddle point 
of the energy, where the magnetization is perpendicular to the easy axis. The 
magnetization vectors of all nodes are shown. Due to the small structure size, 
all spin vectors are almost perfectly parallel to each other, indicating that we 
should expect a good agreement with the analytic result. 

In order to compare the numerical results directly with the analytical results, 
the attempt frequency is plotted as a function of the external field strength as 
shown in [l] (A) . The simulations agree very well with the analytical solution. 
It should be noted, that both the analytical solution and the numerical simula- 
tion are not valid for zero external field, since the expressions for field lowered 
symmetry are used. As a consequence the spike-like increase at zero external 
field is an artifact and has no physical meaning. 

3.2 Elongated particle 

As a second example we calculate the attempt frequency of an elongated par- 
ticle. The simulated geometry of the elongated grain is: length 25 nm, width 
1.6 X 1.6 nm. This shape is similar to the magnetic grains used as the recording 
media in modern hard disks, only slimmer (die lateral dimension of a grain of a 
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hard disc is about 6x6 nm). According to the analytical theory of Ref. |16|the 
cross section area does not influence the attempt frequency. Hence, in order to 
save computational time we model the grain of a recording media with a smaller 
cross section but we can expect that the results can be used for realistic grains. 

The following parameters were used for the graded media grain example: 
damping constant a = 0.2 , anisotropy constant increasing quadratically from 
zero to Ki = 3.6MJ/m^ , magnetic saturation polarization Js = 0.5T and 
exchange constant A = lOpJ/m. The material parameters are the same as in 



directions of the easy-axis and the applied perpendicular external magnetic field 
are shown. The magnetic state at the minimum of the energy is shown on the 
left side of the figure, and the saddle point is shown on the right side. 

As in the previous section the attempt frequency was calculated as a function 
of a perpendicular applied field. In order to verify reliability, the simulations 
were repeated with different finite element meshes. 

Fig. [3] shows the attempt frequency of four different finite element meshes, 
with different numbers of volume elements. As can be seen, the results of the 
attempt frequency do not converge to a single value as the finite element mesh 
size is decreased. Furthermore the obtained value of the attempt frequency is 
one order of magnitude larger than the value reported in Ref [^T]. The reason 
for these discrepancies are discussed in the next section. 

4 Improved method using analytical fits 

The reason for the numerical errors of the attempt frequency can be found in the 
calculation of fioi which consists of (the square root of) a ratio of determinants 
as given by Eq |15| The calculation of the determinant of the symmetric matrix 
H equals the product of all its eigenvalues Ai, . . . , X2N, where N is the 
number of finite element node points. 
Hence, we can write VIq in the form: 



All eigenvalues of H„iin at the minimum are positive. At the saddle point 
there is only one negative eigenvalue. The number of eigenvalues increases 
with increasing number of finite elements. Hence, small numerical errors in 
the calculation of the eigenvalues multiply together to form a considerable total 
error. 

Fig. ^ shows the ratio of the eigenvalues. The x-axis represents the index i 
of the eigenvalues, and the y-axis represents the natural logarithm of the ratios 
of eigenvalues, Xmin,i/ Xsp,i- In Fig- |4]the eigenvalues are sorted according to 
modulus, where the eigenvalue with the index i = 1 has the smallest modulus 
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and i — 2N is the eigenvalue with the largest modulus. As shown in Fig. |4] 
only ratios of eigenvalues with small i are not close to or equal to 1. All the 
ratios of the eigenvalue with i » 100 are very close to one, which results in 
ln{Xmin,i/Kp,i) « for i >> 100. 

The numerical evidence in Fig. |4] illustrate the insensitivity of the saddle 
point to the majority of eigenmodes, and the thermal significance of modes 
with the smallest eigenvalues. As noted in Ref [T5], a domain wall pinned at an 
interface supports a broad range of travelling spin waves that are essentially un- 
perturbed by the wall structure, and a set of modes localized to the domain wall. 
The modes localized to the wall have the smallest eigenvalues, and represent the 
most relevant fluctuations for thermal depinning of the wall. 

Braun has shown that for elongated particles, the logarithm of the ratio of 
eigenvalues scales with (see formula (4.12) of Ref. [T]). 

In order to reduce the numerical error which is introduced by the higher order 
eigenfrequency, we fit the values of lT^{Xmin,i/^sp,i) to the following function. 



/W = «-(^Y^j (16) 

, where a and b are parameters determined by a mean square fit to the 
numerical data for 5 < i < 60. In Fig. [5] the black line represents the ratios 
of eigenvalues, calculated by the numerical simulation, and the red line is the 
fitted function, which values of parameters a and b are given in the inset. The 



calculation of Qq is done by calculating fig = 



-pre . YfN „fit(i) 



We found for a wide range of c (10 < c < 100) that the value of fip is not 
significantly infiuenced by the actual choosen value of the paramter c. 

The calculation of the attempt frequency of the graded media grain is re- 
peated with the improved method using the fit function. Fig. [6] shows that the 
mesh dependency of the results is decreased by this method. If the value of the 
attempt frequency is extrapolated for an external field approaching zero we can 
estimate the attempt frequency to be in the order of /o ^ 1850GiJz ± 65QGHz. 
Comparing this estimate of the attempt frequency with the estimate of the at- 
tempt frequency obtained by Langevin dynamic simulations (In Ref. |21| for 
the same system a value of /o ~ 1638GHz ± AQGHz is estimated), despite 
the various assumptions which are used in both methods, a considerable good 
agreement can be found. 

From an application point of view it is interesting to compare the attempt 
frequency of a grain with graded anisotropy to the attempt frequency of grain 
for which anisotropy is uniform throughout. Fig. [7]shows the attempt frequency 
of a single phase grain, where all the material parameters are the same as for 
the graded media grain of Fig. [6] except that the anisotropy is a constant 
Ki ~ 3.6MJ/m^. This value equals the maximum anisotropy constant of the 
graded media grain of the previous results. 
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5 Results and discussion 



A numerical implementation of the transition state theory was described with 
application to two example problems: reversal of a single phase and graded me- 
dia magnetic grains. The implementation makes use of the micromagnetic pack- 
age FEMME. This method allows for the calculation of the long term thermal 
stability of magnetic nanostructures without any free parameters. TIk^ input 
parameters of the model are the exchange constant, the anisotropy constant the 
spontaneous saturation magnetization, the damping constant and the geometry 
of the magnet. The advantage of the presented method over Langevin-dynamic 
simulation is that the thermal stability of large magnetic nanostructures can be 
calculated. Furthermore Langevin-dynamic simulations are restricted to time 
scales of several nano-seconds due to computational effort. 

In previous work it was shown that so called exchange spring media exhibit 
superior writeability as compared to single phase media, without lowering the 
energy barrier [?]. The energy barrier was assumed to determine the thermal 
stability. In this paper we show the exact relation of the energy barrier and the 
thermal stability by calculating the attempt frequency. In the present paper we 
show that the attempt frequency and energy barrier aare comparable for a single 
phase media and a graded media grain with single phase grain anisotropy is the 
same as the rnaxinmm anisotropy in the graded media grain. Hence, it was 
shown that these two media types indeed have very similar thermal stabilities 
but the coercive field for the graded media grain is about a factor of seven 
smaller than the single phase grain. 

The financial support of the FWF projects P20306-N16 and SFB ViCoM 
(F4112-N13) is acknowledged. 
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Figure 1: (color online) (A) the numerically calculated attempt frequency of 
a small magnetic particle is compared with the analytical result as a function 
of the external applied field (H field applied along y-axis, easy axis parallel to 
X-axis). (B) The magnetic state at the minimum (left) and at the saddle point 
(right) is shown. 
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Figure 2: (color online) model of the graded media grain which is soft magnetic 
at the top and hard magnetic at the bottom, (left) energy minimum state (right) 
saddle point configuration. 
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Figure 3: (color online) The attempt frequency as a function of the external 
field strength is shown for a graded media grain. The external field is applied 
perpendicular to the easy axis. The results are shown as a function of the 
number of finite elements. 



2.5 p 
2.0 
1.5 
t 1.0 
0.5 
0.0 - 



°'^0 100 200 300 400 500 600 700 800 

i {index of A) 

Figure 4: Logarithmic ratios of eigenvalues as a function of the eigenvalue index 
i. (number of finite elements= 1316 , Hext is 0.02 T.) 
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Fitting function is: 7.77076274264 ■(;/(! )) 



2 54267558738 
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Figure 5: Contributions to calculation of Hq: up to index 50 the original data is 
taken (black line), and from index 50 onwards the values of the fitting function 
are used (red line). 
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Figure 6: Same as Fig. [3] except that a fitting function from eigenvalue index 
30 onwards is used for the calculation of Hq. 
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Figure 7: Attempt frequency of a single phase media calculated for two different 
mesh sizes (For the calculation of VIq a fitting function from eigenvalue index 30 
onwards is used). 
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